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The heat capacity of iron isotopes is calculated within the interacting shell 
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Q>^ ' model using the complete {pf + 0^9/2 )-shell. We identify a signature of the 

^D ' pairing transition in the heat capacity that is correlated with the suppression 

of the number of spin-zero neutron pairs as the temperature increases. Our 

Q ' results are obtained by a novel method that significantly reduces the statis- 

(Ih ' tical errors in the heat capacity calculated by the shell model Monte Carlo 



approach. The Monte Carlo results are compared with finite-temperature 



j^ ■ Fermi gas and BCS calculations. 
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Pairing effects in finite nuclei are well known; examples include the energy gap in the 
spectra of even-even nuclei and an odd-even effect observed in nuclear masses. However, less 
is known about the thermal signatures of the pairing interaction in nuclei. In a macroscopic 
conductor, pairing leads to a phase transition from a normal metal to a superconductor below 
a certain critical temperature, and in the BCS theory [|I| the heat capacity is characterized 
by a finite discontinuity at the transition temperature. As the linear dimension of the 
system decreases below the pair coherence length, fluctuations in the order parameter become 
important and lead to a smooth transition. The effects of both static fluctuations p|,|^ 
and small quantal fluctuations IQ have been explored in studies of small metallic grains. 
A pronounced peak in the heat capacity is observed for a large number of electrons, but 
for less than ~ 100 electrons the peak in the heat capacity all but disappears. In the 
nucleus, the pair coherence length is always much larger than the nuclear radius, and large 
fluctuations are expected to suppress any singularity in the heat capacity. An interesting 
question is whether any signature of the pairing transition still exists in the heat capacity 
of the nucleus despite the large fluctuations. When only static and small- amplitude quantal 
fluctuations are taken into account, a shallow 'kink' could still be seen in the heat capacity 
of an even-even nucleus [^. This calculation, however, was limited to a schematic pairing 
model. Canonical heat capacities were recently extracted from level density measurements 
in rare-earth nuclei and were found to have an S-shape that is interpreted to represent 
the suppression of pairing correlations with increasing temperature. 

The calculation of the heat capacity of the finite interacting nuclear system beyond the 
mean-field and static-path approximations is a difficult problem. Correlation effects due to 
residual interactions can be accounted for in the framework of the interacting nuclear shell 
model. However, at finite temperature a large number of excited states contribute to the heat 
capacity and very large model spaces are necessary to obtain reliable results. The shell model 
Monte Carlo (SMMC) method [^^ enables zero- and finite-temperature calculations in large 
spaces. In particular, the thermal energy E{T) can be computed versus temperature T and 
the heat capacity can be obtained by taking a numerical derivative C = dE/dT. However, 



the finite statistical errors in E{T) lead to large statistical errors in the heat capacity at 
low temperatures (even for good-sign interactions). Such large errors occur already around 
the pairing transition temperature and thus no definite signatures of the pairing transition 
could be identified. Furthermore, the large errors often lead to spurious structure in the 
calculated heat capacity. Presumably, a more accurate heat capacity can be obtained by 
a direct calculation of the variance of the Hamiltonian, but in SMMC such a calculation 
is impractical since it involves a four-body operator. The variance of the Hamiltonian has 
been calculated using a different Monte Carlo algorithm f^, but that method is presently 
limited to a schematic pairing interaction. Here we report a novel method for calculating 
the heat capacity within SMMC that takes into account correlated errors and leads to much 
smaller statistical errors. Using this method we are able to identify a signature of the 
pairing transition in realistic calculations of the heat capacity of finite nuclei. The signature 
is well correlated with the suppression in the number of spin-zero pairs across the transition 
temperature. 

The Monte Carlo approach is based on the Hubbard- Stratonovich (HS) representation 
of the many-body imaginary-time propagator, e~^^ = J D[a]GfjUa, where f3 is the in- 
verse temperature, G^ is a Gaussian weight and U„ is a one-body propagator that de- 
scribes non-interacting nucleons moving in fluctuating time- dependent auxiliary fields a. 
The canonical thermal expectation value of an observable O can be written as (O) = 
/D[fT]Go-Tr (OUo-)/ J D[a]GrjTrU„, where Tr denotes a canonical trace for N neutrons and 
Z protons. We can rewrite 

^^^ ^ ([Tr(0^.)/Tr[/.]«l>.)^ ^ ^^^ 

where $0- = TtU(j/\TtU„\ is the Monte Carlo sign, and we have used the notation {X^)w = 
J D[a]W{a)Xj J D[a]W{a) with W{a) = G^\TtU^\. In SMMC we divide the imaginary- 
time interval {0,(3) into Nt time slices of length A/5 = P/Nt, and sample the fields cr(r„) 
fields at all Nt time slices r„ = nA/? according to the positive-definite weight function W{a). 
Each quantity {X^)^r in Eq. (|lD is then estimated as an arithmetic average over the chosen 



samples. 

In particular the thermal energy can be calculated as a thermal average of the Hamilto- 
nian H . The heat capacity C = —(3'^dE/d(3 is then calculated by estimating the derivative 
as a finite difference 

c = -^^^(^^'^\-J(^-'^Kom^. (2) 

At low temperatures, E{(3) decreases slowly with (3 and even small errors in E{(3) lead to 
relatively large statistical errors in C. Conventionally, the calculation of E{f3) for each (3 is 
done by a new Monte Carlo sampling of the a fields and consequently the energies E{f3 — 6f3) 
and E{i3 + 6(3) in (|^) are uncorrelated. However, if the calculation of both E{j3 ± 6(3) can 
be done using a common set of sampling fields, the correlated errors of C are expected to 
be smaller. 

The energies E{(3 ± 6(3) are calculated from 

p.. I e.^ JD[a'^]GMP±S(3)Tv[HUMP±5(3)] 

^'^ ^^ JD[a^]G^±{(3±6(3)TTUMf^±6(3) ' ^^ 

where the corresponding a fields are denoted by (j"^. To have the same number of time 
slices Nt in the discretized version of (0) as in the original HS representation of E{(3), we 
define modified time slices A(3± by NtA(3± = /3 ± 6(3. We next change integration variables 
in (0) from o"^ to a according to o"^ = {AP/AP±)^''^a, so that the Gaussian weight is left 



unchanged G^± {(3 ± 6(3) = exp 



an 



Y.\W\{at{rrO)^AP± =exp - E ||t^a|K(T:„))'A/? 



an 



Gcr{(3) {va are the interaction 'eigenvalues', obtained by writing the interaction in a quadratic 
form J^a^aPaf^, whcrc Pa are one-body densities). Rewriting (^ using the measure D[a] 
(the Jacobian resulting from the change in integration variables is constant and canceled 
between the numerator and denominator), we find 

Tr HU^± (/3±5/3) Tr U^± (/3±5/3) ^ 



Ei(3±6(3) = / T.a,^(/3±./.) . V = ^- (4) 

The heat capacity in (|^) is calculated from G = —(3'^{26(3)^^{H^/Z^ — H^jZJ). Since the 
same set of fields a is used in the calculation of both E{(3^6(3), we expect strong correlations 



among the quantities H± and Z±, which would lead to a smaller error for C. The covariances 
among H^ and Z± as well as their variances can be calculated in the Monte Carlo and used 
to estimate the correlated error of the heat capacity. 

We have calculated the heat capacity for the iron isotopes 52-62pg ^gj^g i^j-^g complete 
{pf + 0(?9/2)-shell and the good-sign interaction of Ref. [^. Fig. |I| demonstrates the sig- 



nificant improvement in the statistical Monte Carlo errors. On the left panel of this figure 
we show the heat capacity of ^^Fe calculated in the conventional method, while the right 
panel shows the results from the new method. The statistical errors for T ~ 0.5 — 1 MeV 
are reduced by almost an order of magnitude. The results obtained in the conventional 
calculation seem to indicate a shallow peak in the heat capacity around T ~ 1.25 MeV, but 
the calculation using the improved method shows no such structure. 

The heat capacities of four iron isotopes ^^~^'^Fe, calculated with the new method, are 
shown in the top panel of Fig. Q. The heat capacities of the two even-mass iron isotopes 
(^^Fe and ^^Fe) show a different behavior around T ~ 0.7 — 0.8 MeV as compared with the 
two odd- mass isotopes (^^Fe and ^"^Fe). While the heat capacity of the odd- mass isotopes 
increases smoothly as a function of temperature, the heat capacity of the even-mass isotopes 
is enhanced for T ~ 0.6 — 1 MeV, increasing sharply and then leveling off, displaying a 
'shoulder.' This 'shoulder' is more pronounced for the isotope with more neutrons (^^Fe). 
To correlate this behavior of the heat capacity with a pairing transition, we calculated the 
number of J = nucleon pairs in these nuclei. A J = pair operator is defined as usual by 

^ = Z^ ~J===Z^jania"'ja-ma ' \^) 

a,m„>0 dja + 1/2 



where ja is the spin and ma is the spin projection of a single-particle orbit a. Pair-creation 
operators of the form (||) can be defined for protons (Aj^^), neutrons (A]^„), and proton- 
neutrons (A],„). The average number (A^A) of J = pairs (of each type) can be calculated 
exactly in SMMC as a function of temperature. The bottom panel of Fig. |^ shows the 
number of neutron pairs (Ajj„A„„) for ^^~^^Fe. At low temperature the number of neutron 
pairs for isotopes with an even number of neutrons is significantly larger than that for 



isotopes with an odd number of neutrons. Furthermore, for the even-mass isotopes we 
observe a rapid suppression of the number of neutron pairs that correlates with the 'shoulder' 
observed in the heat capacity. The different qualitative behavior in the number of neutron 
pairs versus temperature between odd- and even-mass iron isotopes provides a clue to the 
difference in their heat capacities. A transition from a pair-correlated ground state to a 
normal state at higher temperatures requires additional energy for breaking of neutron pairs, 
hence the steeper increase observed in the heat capacity of the even-mass iron isotopes. 
Once the pairs are broken, less energy is required to increase the temperature, and the heat 
capacity shows only a moderate increase. 

It is instructive to compare the SMMC heat capacity with a Fermi gas and BCS calcula- 
tions. The heat capacity can be calculated from the entropy using the relation C = TdS/dT. 
The entropy S of uncorrelated fermions is given by 

S{T) = -E[/aln/a + (1 - /a)ln(l - /.)] , (6) 

a 

with fa being the finite-temperature occupation numbers of the single-particle orbits a. 
Above the pairing transition-temperature T^, fa are just the Fermi-Dirac occupancies fa = 
[1 + e'^*^'^""'^)]^^, where fx is the chemical potential determined from the total number of 
particles and €„ are the single-particle energies. Below Tc, it is necessary to take into account 
the BCS solution which has lower free energy. Since condensed pairs do not contribute to 
the entropy, the latter is still given by (^) but fa are now the quasi-particle occupancies [|l| , 

-E'a = Y (ca — yu)^ + A^ are the quasi-particle energies, where the gap A(T) and the chemical 
potential fi{T) are determined from the finite-temperature BCS equations. In practice, we 
treat protons and neutrons separately. 

We applied the Fermi gas and BCS approximations to estimate the heat capacities of the 
iron isotopes. To take into account effects of a quadrupole-quadrupole interaction, we used 
an axially deformed Woods-Saxon potential to extract the single-particle spectrum €„ |]ll 



A deformation parameter S for the even iron isotopes can be extracted from experimental 
B{E2) values. However, since B{E2) values are not available for all of these isotopes, we 
used an alternate procedure. The excitation energy Ex{2f) of the first excited 2+ state 
in even-even nuclei can be extracted in SMMC by calculating (J^)/3 at low temperatures 
and using a two-state model (the 0"*" ground state and the first excited 2"*" state) where 
( J2)^ ^ 6/(1 + e^^^^^i V5) IH]. The excitation energy of the 2i state is then used in the 
empirical formula of Bohr and Mottelson H r^ = (5.94 ± 2.43) x 10^^E-^{2t)Z-^A^/^ to 
estimate the mean 7-ray lifetime r^ and the corresponding B{E2). The deformation pa- 
rameter 5 is then estimated from B(E2) = [(3/47r)ZerQy4^/^5]^/5. We find (using ro = 1.27 
fm) 6 = 0.225, 0.215, 0.244, 0.222, 0.230, and 0.220 for the even iron isotopes ^^Pe - ^^Fe, 



respectively. For the odd- mass iron isotopes we adapt the deformations in Ref. [Q. The 
zero-temperature pairing gap A is extracted from experimental odd-even mass differences 
and used to determine the pairing strength G (needed for the finite temperature BCS solu- 
tion) . 

The top panels of Fig. ^ show the Fermi-gas heat capacity (dotted-dashed lines) for 
^^Fe (right) and ^°Fe (left) in comparison with the SMMC results (symbols). The SMMC 
heat capacity in the even- mass ^°Fe is below the Fermi-gas estimate for T ^ 0.5 MeV, 
but is enhanced above the Fermi gas heat capacity in the region 0.5 ^ T ^ 0.9 MeV. 
The line shape of the heat capacity is similar to the S'-shape found experimentally in the 
heat capacity of rare-earth nuclei 0. We remark that the saturation of the SMMC heat 
capacity above ~ 1.5 MeV (and eventually its decrease with T) is an artifact of the finite 
model space. The solid line shown for ^°Fe is the result of the BCS calculation. There are 
two 'peaks' in the heat capacity corresponding to separate pairing transitions for neutrons 
{T^ ^ 0.9 MeV) and protons (T^ ^ 1.2 MeV). The finite discontinuities in the BCS heat 
capacity are shown by the dotted lines. The pairing solution describes well the SMMC 
results for T < 0.6 MeV. However, the BCS peak in the heat capacity is strongly suppressed 
around the transition temperature. This is expected in the finite nuclear system because 
of the strong fluctuations in the vicinity of the pairing transition (not accounted for in the 



mean-field approach). Despite the large fluctuations, a 'shoulder' still remains around the 
neutron-pairing transition temperature. 

The bottom panels of Fig. ^ show the number of spin-zero pairs versus temperature in 
SMMC The number oi p-p and n-p pairs are similar in the even and odd-mass iron isotopes. 
However, the number of n-n pairs at low T differs signiflcantly between the two isotopes. 
The n-n pair number of ^''Fe decreases rapidly as a function of T, while that of ^^Fe decreases 
slowly. The S-shape or shoulder seen in the SMMC heat capacity of ®°Fe correlates well with 
the suppression of neutron pairs. 

Fig. H shows the complete systematics of the heat capacity for the iron isotopes in 
the mass range A = 52 — 62 for both even- mass (left panel) and odd- mass (right panel). 
At low temperatures the heat capacity approaches zero, as expected. When T is high, the 
heat capacity for all isotopes converges to approximately the same value. In the intermediate 
temperature region (T ~ 0.7 MeV), the heat capacity increases with mass due to the increase 
of the density of states with mass. Pairing leads to an odd-even staggering effect in the mass 
dependence (see also in Fig. ^ where the heat capacity of an odd-mass nucleus is signiflcantly 
lower than that of the adjacent even-mass nuclei. For example, the heat capacity of ^^Fe is 
below that of both ^^Fe and ^^Fe. The heat capacities of the even-mass ^^Fe, ^''Fe, and ^^Fe 
all display a peak around T ~ 0.7 MeV, which becomes more pronounced with an increasing 
number of neutrons. 

In conclusion, we have introduced a new method for calculating the heat capacity in 
which the statistical errors are strongly reduced. A systematic study in several iron isotopes 
reveals signatures of the pairing transition in the heat capacity of flnite nuclei despite the 
large fluctuations. 
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FIG. 1. The SMMC heat capacity of ^^Fe. The left panel is the result of conventional SMMC 
calculations. The right panel is calculated using the improved method (based on the representation 
@ where a correlated error can be accounted for). 
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FIG. 2. Top panel: the SMMC heat capacity vs. temperature T for ^^Fe (open circles), ^^Fe 
(solid diamonds), ^^Fe (open squares), and ^^Fe (solid triangles). Bottom panel: the number of 
J = neutron pairs versus temperature for the same nuclei. 
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FIG. 3. Top: Heat capacity versus T for ^"Fe (left) and ^^Fe (right). The Monte Carlo results 
are shown by symbols. The dotted-dashed lines are the Fermi gas calculations, and the solid line 
(left panel only) is the BCS result. The discontinuities (dashed lines) correspond to a neutron 
{Tc ^ 0.9 MeV) and proton (Tc ~ 1.2 MeV) pairing transition. Above the pairing-transition 
temperature, the BCS results coincide with the Fermi gas results. Bottom panels: The number of 
J = n-n (circles), p-p(squares), and n-p(diamonds) pairs vs. T for ^'^Fe (left) and ^^Fe (right). 
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FIG. 4. The heat capacity of even-even (left panel) and odd-even (right panel) iron isotopes. 
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